PROGRAM traj
IMPLICIT NONE
REAL :: v,theta,k,m
REAL :: x, y, t, total_time
REAL :: vx, vy, vx0, vy0, x0, y0, t0
REAL,PARAMETER :: PI=3.14159265359, g=9.8
INTEGER :: eqn
WRITE (*,*), 'Enter velocity:'; READ *, v  !Reading the values from the user
WRITE (*,*), 'Enter Angle of launch from horizontal (in degree):'; READ *, theta
WRITE (*,*), 'Enter Mass (in kg):'; READ *, m

theta=theta*PI/180
total_time=v*SIN(theta)/g

OPEN(1,FILE='q5_drag0.dat')
k=0
vx0=v*COS(theta)
vy0=v*SIN(theta)
x0=0
y0=0

DO t=0.0,2*total_time,0.5
  x=vx0*t
  y=vy0*t - 0.5*g*t*t
  WRITE (1,*),x,y
END DO
CLOSE(1)

OPEN(2,FILE='q5_drag0.1.dat')
k=0.1
vx0=v*COS(theta)
vy0=v*SIN(theta)
x0=0
y0=0
DO t=0.0,total_time,0.5
  x=vx0*m/k - (vx0*m/k)*EXP(-1*k/m*t)
  y=(vy0+m*g/k)*m/k*(1-EXP(-1*k*t/m)) - m*g*t/k
  vy=(vy0+m*g/k)*EXP(-1*k*t/m) - m*g/k
  IF (vy>=0) THEN
    WRITE (2,*),x,y
  ELSE
   EXIT
  ENDIF
END DO
t0=t
y0=y
PRINT *,'t=',t0,total_time
DO t=0.0,total_time,0.5
  x=vx0*m/k - (vx0*m/k)*EXP(-1*k/m*(t+t0))
  y=g/k*(t+m/k*(1-EXP(k*t/m)))
  y=y0+y
  IF(y>0) THEN
    WRITE (2,*),x,y
  ELSE
    EXIT
  ENDIF
END DO
CLOSE(2)

OPEN(2,FILE='q5_drag0.2.dat')
k=0.2
vx0=v*COS(theta)
vy0=v*SIN(theta)
x0=0
y0=0
DO t=0.0,total_time,0.5
  x=vx0*m/k - (vx0*m/k)*EXP(-1*k/m*t)
  y=(vy0+m*g/k)*m/k*(1-EXP(-1*k*t/m)) - m*g*t/k
  vy=(vy0+m*g/k)*EXP(-1*k*t/m) - m*g/k
  IF (vy>=0) THEN
    WRITE (2,*),x,y
  ELSE
   EXIT
  ENDIF
END DO
t0=t
y0=y
PRINT *,'t=',t0,total_time
DO t=0.0,total_time,0.5
  x=vx0*m/k - (vx0*m/k)*EXP(-1*k/m*(t+t0))
  y=g/k*(t+m/k*(1-EXP(k*t/m)))
  y=y0+y
  IF(y>0) THEN
    WRITE (2,*),x,y
  ELSE
    EXIT
  ENDIF
END DO
CLOSE(2)

PRINT *,'WROTE DATA TO FILES'
END PROGRAM

